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Abstract.  The  paper  addresses  the  a-posteriori  error  estimates  of  modelling 
heat  convection  problems  in  a  plate.  It  derives  estimators  which  are 
essentially  local  (of  the  size  of  the  thickness  of  the  plate).  The  estimator 
is  guaranteed  to  be  an  upper  bound.  The  lower  bound  of  the  error  is  also 
given.  The  asymptotic  exactness  of  the  estimator  with  respect  to  the 
thickness  of  the  plate  and  with  respect  to  the  order  of  the  model  is  proven. 
The  adaptive  procedure  based  on  this  a-posterlori  error  estimator  is 
proposed.  Numerical  examples  are  given. 
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1.  Introduction 


During  recent  years  a  large  progress  has  been  achieved  in  the  field  of 
the  a  posteriori  error  estimations  and  adaptive  procedures  in  the  finite 
element  method.  The  basic  ideas  have  been  introduced  in  [1-3]  and  since  then 
many  error  estimators  appeared  in  the  literature. 

The  basis  of  these  estimators  are  essentially  of  the  following  types. 

a)  Residual  estimators  which  employ  the  residuals  of  the  finite 
element  solution.  For  more  about  these  we  refer  among  others  to  [1-13]. 

These  estimators  often  are  utilizing  the  complementary  energy  principles 
(see  e.g.  [4,5,11,13]). 

b)  Flux-projection  operators.  These  estimators  are  based  on  the 
idea  of  smoothing  the  fluxes  (which  are  discontinuous)  and 

comparison  with  the  fluxes  computed  from  the  original  finite  element  solution. 
For  more  we  refer  to  [14-18]. 

c)  Other  categories  include  extrapolation  estimators  (see  e.g.  [19]) 
and  interpolation-error  bounds,  see  e.g.  [20-21]. 

In  [22-24]  the  comparison  of  the  performance  of  various  major  error 
estimators  and  their  robustness  is  given  (see  also  there  for  an  extensive 
list  of  references). 

The  a  posteriori  error  estimators  are  the  basis  of  the  adaptive 
procedures.  The  indicators  (which  are  the  error  estimators  in  the  single 
elements)  govern  then  the  design  of  the  mesh. 

We  underline  that  the  error  is  understood  here  as  the  difference 
between  the  exact  solution  of  the  differential  equation  under  consideration 
and  its  finite  element  approximation.  The  addressed  differential  equation 
describes  the  particular  model,  say,  of  the  plate  and  the  computed  error 
characterizes  the  accuracy  of  the  finite  element  method  with  respect  to 
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the  exact  solution  of  the  selected  plate  model  and  not  with  respect  to  the 
plate  itself.  Nevertheless  the  (exact)  solution  of  the  model  (plate  or 
shell)  problem  has  to  be  understood  as  the  approximate  solution  of  the 
3~dimensional  formulation. 

Many  models  were  proposed  in  the  literature.  Ve  refer  for  example  to 
the  survey  [25].  It  is  shown  in  [26]  that  various  plate  models  could  lead  to 
significantly  different  results  especially  in  the  neighborhood  of  the  boundary 
and  in  the  presence  of  unsmooth  input  data.  Hence  the  (exact)  solution  of  the 
model  problem  has  to  be  understood  as  an  approximation  of  the  (exact)  solution 
of  another  "higher"  problem. 

In  the  example  we  mentioned,  the  "higher"  problem  is  the  three 
dimensional  model. 

Similarly  as  in  any  approximate  method  a  hierarchy  of  models  has  to  be 
available  and  be  such  that  it  allows  to  obtain  the  solution  of  the  original 
"higher"  model  with  a  priori  given  accuracy.  Practically  this  means  that  an 
a-posteriorl  error  estimator  of  the  error  of  the  model  in  comparison  with  the 
exact  solution  of  the  higher  model  should  exist  and  an  adaptive  procedure  for 
the  optimal  model  selection  should  be  available. 

In  practice  the  exact  solution  of  the  model  is  impossible  to  find 
either.  Hence  we  have  to  combine  both  errors — the  error  of  the  model  and  of 
the  finite  element  solution.  In  the  adaptive  procedure  both  errors  have  to 
be  controlled. 

In  contrast  to  the  a  posteriori  error  estimation  techniques  and  adaptive 
approaches  of  the  finite  element  method  as  described  briefly  earlier,  the  a 
posteriori  error  analysis  of  the  models,  and  the  problem  of  adaptive  modelling 
has  not  been  addressed  in  the  literature  until  now.  In  this  paper  we  address 
this  problem  for  the  solution  of  a  model  heat  conduction  problem  in  laminated 
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plates.  Let  us  mention  that  this  problem  has  also  other  physical 
interpretations.  The  theoretical  results  presented  here  are  based  on  [27], 

We  concentrate  in  this  paper  on  the  error  of  the  models,  and  its 
adaptive  selection.  We  will  only  consider  the  case  when  the  exact  solution 
of  the  models  is  available.  Hence  we  will  assume  that  the  error  of  the 
finite  element  solution  is  negligible  in  comparison  with  the  model  error. 

The  paper  is  organized  as  follows.  After  introducing  some  notation  and 
the  three-dimensional  boundary  value  problem  in  Section  2,  we  put  in  Section 
3  the  hierarchical  modelling  in  perspective  with  other  approaches  to  derive 
reduced  models.  In  Section  4  we  collect  basic  properties  of  the  hierarchical 
models  from  [27].  Section  5  contains  the  derivation  of  the  a  posteriori 
estimator  for  the  modelling  error  and  demonstrates  its  asymptotic  exactness 
as  d — >0,  as  well  as  its  spectral  exactness  as  q— »».  In  Section  6  we 
illustrate  our  theory  by  a  simple  numerical  example,  demonstrating  in 
particular  that  our  theoretical  bounds  of  the  effectivity  indices  are  sharp. 
While  Section  6  dealt  with  a  smooth  solution,  we  investigate  in  Section  7  the 
performance  of  our  estimator  in  the  presence  of  boundary  layers.  In  Section  8 
we  address  the  adaptive  selection  of  the  model. 
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R+  =  {(x1,x2,x3) j(x1(x2)  €  u,  X3  =  ±  d/2>. 

Often  we  will  write  x  =  (x^Xg),  x  =  (Xj,x2>x3)  =  (x,x3). 
As  a  model  problem  we  will  consider  the  problem 


A  u 

=  0 

on 

n, 

(2.  1) 

u 

=  0 

on 

r, 

D  u 
n 

=  f 

on 

R±, 

where  D  is  the  exterior  unit  normal  derivative, 
n 

Let 

H  ■  (u  e  HX(n)  |  u  =  0  on  D 

and  define  the  bilinear  form  B{  • ,  • )  :  HxH — >!R  and  the  functional 
F(  ■  ) :  H  — »IR 

(2.2)  B(u.u)  =  Jvu-Vudx. 

a 

(2.3)  F(u)  =  Jf (x1,x2)(u(x1,x2>|)  +  u(x1,x2>  -|))dx. 

n 

Then  the  weak  form  of  (2.1)  reads: 

Find  u  e  H  such  that 

(2.4)  B(u, v)  =  F(u) ,  V  u  e  H. 

There  exists  a  unique  solution  of  the  problem  (2.2)-(2.4)  provided  that 
f  €  L-2 ( c*>)  (the  assumption  on  f  can  be  weakened).  Further  we  define  the 
energy  norm 

(2.5)  |»  u  III2  =  B(u,u)2. 

We  remark  that  the  problem  (2.1 )— (2. 3)  respectively  (2.4)  describes  the 
heat  conduction  problem  in  the  plate  n.  We  assumed  that  the  material  of  the 
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plate  is  homogeneous;  nevertheless  the  theory  Is  valid  with  small  charges 
for  the  laminated  plates,  too  (see  [27]). 


3.  The  modelling  problem 

The  problem  of  the  modelling  means  here  a  reduction  of  the  three 
dimensional  problem  (1.1)  to  a  two  dimensional  one.  There  are  three  basic 
types  of  such  reductions  (which  are  available  not  only  for  the  heat 
conduction  problem  but  also  for  plates  and  shells).  They  are 

1)  Physical  derivation, 

2)  Asymptotic  derivation, 

3)  Hierarchical  dimensional  reduction. 

Let  us  describe  the  main  ideas  of  these  three  approaches. 

1)  Physical  derivation 

The  problem  (2.1)  is  first  written  in  the  form 

(3.  la)  Vu  =  <r  in  fl 

(3.  lb)  div  <r  =  0  in  Q 

(3.1c)  T-n  =  r  on  R± 

(3. Id)  u  =  0  on  T. 

Here  <r  defined  in  (3. lal  has  the  meaning  of  fluxes,  (3. lb)  and  (3. lc)  have 
the  meaning  of  heat  balance.  For  small  d  we  assume  that  the  flux  <r^  is 
negligible  so  that  the  temperature  u  through  the  thickness  is  constant. 
Then  the  heat  balance  through  the  thickness  reads 

d  div  ( <r1 ,  «r2 )  =  -2  f . 

This  leads  to 

(3.2a)  -A.u  =  in  n, 

x  d 

(3. 2b)  u  =  0  on  ?. 
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The  second  phase  Is  to  derive  an  Improved  model.  We  will  follow  the 
reasoning  given  in  [29]  (when  applied  to  our  problem). 


Assume  that  the  solution  u  of  (3.2)  is  available  and  denote  it  by  Uj. 
Then  the  corresponding  flux  is 


(3.3) 


(<rr^20))  =  7V 


The  heat  balance  equation  in  11  now  reads 


(3.4) 


div.  <r  +  =  0 

x  dx3 


which  suggests 


(3.5) 


a  o-3  .. 

5-5  =  -div.  <r 
3x3  x 


(3.6) 


<r_  =  -x„  div.  <r  =  -x.A.u  =  2 -rX_. 
J  J  x  3  x  a  J 


<r3(x)  satisfies  then  the  originally  prescribed  boundary  condition 


(3  7) 


°3(xl ' x2’  ±d)  “  f 


and  the  improved  flux  is  now  (o-j ,  <r2>  o-^) .  It  satisfies  the  balance  through 
the  thickness  and  the  boundary  condition  at  x^  =  ±d/2.  Now  using 


(3.8) 


we  get  the  improved  solution  u 


Sx3  =  °3 


(3.9) 


~  12 
u(xltx2,x3)  =  Jf(xl’x2)x3  +  ■ 


Taking  the  heat  balance  through  the  thickness  we  get 

d  Au  =  d3A.f  +2f  +  d  A.^  =  0 

12  x  x 


and  hence 


(3. 10) 


-d  A_t//  =  2f  +  d3A.f  in  u. 

x  12  x 
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The  boundary  condition  at  y  has  to  be  understood  In  a  weak  sense 


J  u(Xj ,  x^,  Xg)  dXg  =  0  for  x^.x^y  which  yields 

^(xltx2)  =  d<Jf(x1,x2),  ( x ^ , x2 )  e  y 

Now  we  can  continue  In  the  same  vein  and  construct  higher  order  approximations 
(in  d)  in  addition  to  u  in  (3.9). 

The  approaches  of  this  type  are  used  very  often  in  the  engineering 
literature.  Nevertheless  their  relations  to  the  original  problems  are  not 
direct.  In  fact,  in  (30]  p.  33  related  to  theories  of  this  type  for  the 
plate  problems  one  finds  the  statement  "the  refined  theories  do  not  yield 
reliable  information  from  the  standpoint  of  3  dimensional  problems. 
Nevertheless,  these  theories  specify  the  principal  stress  of  the  plate". 

2.  Asymptotic  derivation.  For  more  about  asymptotic  analysis  we  refer 
e. g.  to  [31].  We  will  present  here  only  the  main  idea. 

First  the  problem  (2.1)  is  scaled  to  unit  thickness  d  =  1  by  the 
substitution  £  =  (x^,x2),  £3  =  ^2.  Then  the  (2.1)  takes  the  form 

2 

(3. lia)  -A.u  -1-^=0 

?  d2  a?3 

with  the  boundary  condition 

(3.11b)  f0r  53  =  ±? 

(3. 11c)  u  =  0  for  £  e  y. 

We  look  for  a  solution  in  the  form 

(3.12)  u(£,£3)  *  5  u(0)(£)  +  du(1)(£)  +d3u(2)(£)  +  •••  . 

Inserting  (3.12)  into  (3.11)  and  matching  powers  of  d,  we  get  a  sequence  of 
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differential  equations  for  u^^.  The  boundary  condition  (3.11c)  at  T  is 
generally  not  satisfied  by  (3.12)  and  must  be  accounted  for  by  additional 
boundary  layers  of  three-dimensional  character. 

The  main  problem  with  this  and  the  previous  approach  is  that  they  lead 
ultimately  to  the  differentiation  of  the  data  f.  It  is  well  known  that  for 
unsmooth  f  and  fixed  d  an  increase  of  the  number  of  terms  in  (3. 12) 
leads  to  a  decreas  e  of  the  accuracy. 

It  is  also  obvious  that  neither  one  of  the  previous  two  methods  creates 
a  set  (hierarchy)  of  models  by  which  the  exact  solution  of  the  three 
dimensional  problem  can  be  approached  to  any  prescribed  accuracy  at  fixed, 
positive  thickness.  Hence  this  type  of  modelling  is  not  usable  for  our 
purposes. 

3.  The  hierarchic  modelling. 

Denote  by 


(3. 13) 


T  =  {Uj |wj  :  u,  1  <  i  <  n} 


a  collection  of  n  domains  with  piecewise  smooth  boundaries  do.  such  that 

n 

o.  n  o.  =  0  if  i  *  j  and  o  =  U  w  (?  could  be  for  example  a 
i  J  t=l  J 

triangulation  of  w) . 

For  a  vector  of  nonnegative  integers 


(3. 14) 


qn  =  (ql . V-  qi  2  ° 


and  a  sequence  of  linearly  independent  functions  ^  e  H  (-1,1)  we  define 


the  space 

(3. 15) 


qt 

S(?.q):  *  {“  *  H  M„  -  I  U5n<*)*J(??]-  “i  «  T] 

j=0 


Then  S(?,q)  c  H  and  the  (?, q)-mode!  is  the  boundary  value  problem: 
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Find  u(?,q)  €  S(y\q)  such  that 


(3. 16) 


B(u(?,q),v)  *  F(v) ,  Vv  €  S(7\q). 


By  this  we  obtain  a  system  of  elliptic  equations  which  describe  the 
(?,q) -model. 

The  essential  question  is  how  to  select  the  functions  .  This  problem 
was  analyzed  in  (28],  [33].  For  the  laminated  plates  an  analogous  approach 
can  be  used.  In  [28],  [33]  it  has  been  shown  that  a  good  choice  is  (because 
of  the  symmetry  in  x3  of  the  solution) 

(3.17)  ^J(T))  =  P2j(7))’ 

P^(t))  is  the  £-th  Legendre  polynomial.  In  [33]  the  optimal  choices  of  \fi. 
in  dependence  on  various  optimal ’ty  criteria  were  also  analyzed. 


4.  The  basic  properties  of  the  hierarchic  modelling 

Here  we  will  address  only  the  system  (3.17)  of  functions  and  we  will 

2 

deal  here  only  with  the  energy  norm  |j|u8|  =  B(u,u).  In  [27]  we  analyzed  also 

other  norms.  Denoting  by  u  the  exact  solution  and  by  uOP, q)  the  exact 
solution  of  the  (?, q)-model ,  we  define  the  modelling  error 

e(?,q)  =  u-u(f.q) 

and  will  be  interested  in  the  a-posterlori  estimation  of  |||e|||.  In  [27]  the 
following  theorem  has  been  proven. 


Theorem  4.1. 
0  <  d  £  A 

where 


Let  f  €  H^w) 


and  have  compact  support  in  w. 


Then  for 


«e(7\q)ffl  £  C(A)d2N+1/2  ||f|| 

>rN(w) 


and 

u,  but  Is 


N  =  min  qt 

's  the  usual  Sobolev  space.  The  constant 
Independent  of  d  and  f. 


C  depends  on  A  and 
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Let  us  underline  that  the  assumption  about  the  compact  support  of  f  is 

almost  essential.  It  can  be  weakened  but  certain  boundary  compatibility 

conditions  have  to  be  satisfied.  Without  these  conditions  the  solution  u 

has  a  boundary  layer  and  Theorem  4.1  is  not  valid;  the  overall  accuracy  as 

1/2 

d — >0  is  at  best  d  in  general.  Nevertheless,  in  [34]  it  is  shown  that 

the  low  accuracy  is  confined  only  to  a  small  neighborhood  of  y.  In 

2N+1/2 

addition,  the  full  accuracy  of  d  can  be  recovered  by  utilizing  simul¬ 
taneously  different  model  orders  a r  in  the  interior  and  near  y. 

Theorem  4.2.  Assume  that  f  e  H^Cu)  and  0  <  d  <  A.  Let 

where  =  {(x^.x,,)  e  w  |  dist  (Xj.Xg  )  >  t},  w2  =  <‘>S'“l  and  lei  qi  =  N’ 

=  M  i  N.  Then,  for  t  >  C^Nd|6id|  and  M(N,d)  suitable, 

|e(?,q;i  S  C(A)d2N+1/2  ||f || 

Here  depends  only  on  the  Lipschltz  constant  of  y. 

Finally  we  have  Theorem  4.3  which  guarantees  that  we  can  select  a 
sequence  of  models,  the  solution  of  which  converges  to  the  solution  of  the 
three  dimensional  problem  (2.3).  We  have 

Theorem  4.3.  Let  d  >  0  fixed,  f  €  Lg (<»>)■  Then  ||e(?>,q)|||  — >0  an  N  = 

min  q,  — >  ». 
lilSn  1 

The  introduced  hierarchy  is  a  large  family  of  models  which  allows  to 
select  the  models  having  the  following  properties  in  the  major  special  caises 

a)  If  f  is  smooth  and  d  is  small  then  the  particular*  model  leads  to 
an  error  which  is  the  smallest  possible  with  respect  to  d,  i.e.  we  have 
asymptotic  convergence  with  respect  to  d— >0  with  maximal  possible  rate, 
regSLTdless  of  compatibility  conditions  violated  by  the  data. 

b)  In  all  cases,  for  any  error  tolerance  there  exists  a  model  in  the 
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family  which  leads  to  the  error  in  the  prescribed  range. . 

c)  As  will  be  seen  later  an  effective  a-posteriori  error  estimator  of 
the  modelling  error  exists  and  based  on  it  an  adaptive  model  selection 
procedure  can  be  designed. 

The  following  theorem  plays  an  important  role  in  the  theoretical  basis  of 
the  study  of  the  adaptive  hierarchal  modelling. 


Theorem  4.4.  For  every  q)  we  have 


d/2 

J  el 

J-d/2 


efx^.x^.x^)  dXg  =  0  a.e  x  e  u 


where  we  denoted  e(?\q)  =  e. 

We  mention  that  in  fact  we  have  additional  properties  of  e  which  were 
utilized  in  [27]. 

5.  The  a-posterlorl  error  estimation 

Here  we  will  give  some  basic  results  from  [27]  related  to  the 
a-posteriori  estimation  of  the  modelling  error,  when  measured  in  the  energy 
norm.  The  estimates  in  L^-norm  are  also  given  in  [27].  The  a-posteriori 
~rror  estimation  will  be  basis  for  the  adaptive  modelling  .  For  simplicity 
we  will  present  here  the  results  for  f  =  <w}  i.e.  we  assume  that  we  will  not 
partition  w  and  that  q^  =  qQ,  i  -  l,...,n. 

In  the  sections  7  and  8,  we  will  apply  the  results  for  the  general  case. 
Let  0  <  ^(XjXg)  €  W*,so(w)  and 


(5.  1) 


Q:=  max  ||§£-  /2' 

1=1,2  xi' 


L°°(w) 


We  will  be  Interested  in  the  error  measured  in  the  weighted  energy  norm 
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(5.2) 


J  |Ve|2  <p(x1x2) 


dxldx2dx3‘ 


A  typical  choice  is  for  example 


(5.3) 


,0  0.  oc  ,  |  0.  .  (0)  . . 

^(Xj.Xg)  =  exp  —  (|x1~x1l  +  I*2_x2 


where  p  €  (0,1)  and  0  <  a  <  We  can,  of  course,  use  other  choices  of  <p 
too.  Assume  now  that  we  have  computed  the  exact  solution  u(q)  of  the  model 
Then  with 

(5.4)  r(xrx2)  =  f  -  (x1(x2,d/2) 


we  define  the  error  estimator  g  (u(q)) 

q 

1/2 

(5.5)  gq(u(q))  =  [f^rVdXjdXg]  . 

To  assess  the  quality  of  this  estimator  we  will  define  the  effectivity 

index 

g  (u(q) ) 

(5-6)  V'1""1  *  irar  • 

V 

We  will  say  that  g  is  an  upper  (lower)  estimator  if  0  >  1  (0  <  1, 
respectively).  Further  we  will  call  g  proper  with  respect  to  a 

class  T  of  data  f,  if 

(5.7)  0  £  k.  S  0  S  Kn  <  « 

1  q  2 

holds  for  all  data  f  e  T. 

The  estimator  gq  is  asymptotically  resp.  spectrally  exact  on  the  set 

of  the  data  T,  resp.  TQ  if  0  — >  1  as  d— »0+  resp.  0  — >f  as 
u  W  q  q 

q  — »». 
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Let  us  remark  that  as  usual  the  asymptotic  exactness  requires 
considerably  stricter  assumptions  on  the  data  than  the  properness.  Further, 
we  underline  that  the  quality  of  the  estimators  has  to  be  understood  relative 
to  a  class  of  data. 

Let  vis  finally  remark  that  the  weight  can  rapidly  change  in  comparison 
with  d  (see  (5.3)).  This  makes  the  estimator  very  local  and  it  can  be  used 
as  basis  of  the  adaptive  modelling.  In  [27]  we  have  proven  the  following 
theorem 


Theorem  5.1. 


a)  If  f  e  T  =  L2(w)  then 


(5.8) 


0  >  (1-V?  dQ)  (1  +  —  dQ(l  +  —  dQ))  1/2 

it  it 


b)  If  f  e  T0:=  {f  I  M 


0' 


S  |3  <  ®)  ( where  |||r||| 


k  ,<P 


=  f  I  V^r|2^p2dx.dx_)  then 


(5.9) 


0  *  1  +  | 


o2  ♦  q2)11/2. 


(2q+3)  -4  J 

If,  moreover,  <p  si,  then  Q  =  0  and  the  factor  3/2  in  (5.9)  can  be 
replaced  by  1/2. 

Theorem  5.1  shows  that  the  estimator  is  x^,  x ^  proper  with  known  x^  and 

Xg.  If  we  select  the  exponential  weight  <p  given  in  (5.3)  then  the 

estimator  is  asymptotically  exact  on  T0. 

p 

We  see  from  the  theorem  that  we  get  an  upper  estimator  for  a  large  class 
of  data,  namely  we  have  only  to  assume  that  f  e  L^(u).  On  the  other  hand,  to 
get  8-1  we  need  more  smoothness  of  the  solution,  for  example  that  no 
boundary  layer  exists  (this  can  be  achieved  for  example  if  f  is  smooth  and 
has  compact  support  similarly  as  in  Theorem  4.1).  Using  in  (5.8)  and  (5.9) 
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Q  ~  j  we  get  a  practically  local  estimator  of  high  quality. 

Theorem  5.1  concentrates  on  the  case  when  d— >0.  Let  us  now  formulate 
the  theorem  when  d  is  fixed  and  q— In  [27]  we  have  proven 


Theorem  5.2. 


,  1-E 


a)  If  f  e  T_  with  0  =  0(q/d)  ,  e  >  0,  0  independent  of  d,q  then 

P 


C  frfl 

*  l1  +  ci LiJ  J  *  1  =  1-2 


e(q) 


where  Cj  is  independent  of  q  and  d. 

b)  Let  <Pk(xi»x2^  denote  the  k-th  eigenpair  of  -A  in  u>  with 

<p^  |  =  0.  Let  further  Hx^Xg)  =  ^  Pk  ^(Xj.JCg)  and  assume  f  is  such 


that  *  0  and 


k>2 

for  some  e  >  0.  Then 


I  [‘  *  ^ 


IPjjl  -1+e 

iL_  «S  (V2  -  1)  (1  -  q  1  e) 
k 


Ip 


1' 


e(q)  £  (1  -  C2Q)(1  +  2C2Q(1  +  2C2Q)) 


-1/2 


where  Q  is  as  in  (5.1)  and 


C2  = 


and  8  is  the  best  constant  in 


W  <(2q+3)2-4>~E 


Vl1,.,  ,  s  ek- 

L  («) 

In  theorems  5.1  and  5.2  we  analyzed  only  the  error  in  energy  norm. 

2 

Nevertheless,  similar  results  can  also  be  obtained  for  the  L  -norm  estimator 
derived  in  [27]. 


6.  A  simple  example 

The  results  we  have  shown  in  the  previous  sections  hold  also  in  the  case 
that  u  =  (-1,1).  In  the  example  here  we  use  fl  =  w  x  (-d/2,  d/2)  and 
f  =  cos  ^  x^.  Then  obviously 
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(6.  1) 


ch(£  x?) 

u(x  ,x  )  =  cos  (£  X  )  - — — -  . 

In  this  case  we  can  compute  (by  computer  algebra)  the  effectivity  index  0 


explicitly  and  get 

(6.2) 


a2  =  d_  =  ,  ♦  A  .  0(d4). 

"  lei2  “q 


The  coefficients  m  are  listed  in  the  Table  6.1. 

q 


Table  6.1.  The  coefficients  m  in  (6.2) 

q 


q 

0 

1 

2 

3 

4 

5 

6 

m 

q 

240 

360 

936 

1768 

2856 

4200 

5800 

We  see  that  in  fact  the  estimator  is  an  upper  one  and  that  it  is 
asymptotically  exact. 

In  the  Table  6.2  we  report  the  values  of  0^  computed  directly.  The 
finite  element  method  with  very  fine  mesh  and  high  degree  of  elements  p  was 
used  to  ensure  that  we  obtained  an  accurate  solution  of  the  model.  We 
selected  d  =*  1  and  d  =  2,  l.e.  a  very  large  thickness  to  show  the 
effectiveness  of  the  estimator.  In  the  Table  6.7  we  also  report  0 

a2  =  ,  * 

m 

q 

4 

l.e.  we  neglected  the  term  0(d  )  in  (6.2). 

Finally  we  list  the  estimates  based  on  (5.8)  (5.9).  In  our  case  r(x)  = 

Ceos  ^  Xj  and  hence  8  =  In  addition  we  report  the  relative  energy  error 
in  */. 
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Table  6.2.  The  effectivity  indices  and  relative 
errors  for  large  d. 


q 

0(q) 

0(q) 

Upper  b 

Lower  b 

l|e||*/. 

0 

1.0200 

1.0065 

1.32 

1 

40.62 

i 

1.01324 

1.0043 

1.085 

1 

1.057 

2 

1.0052 

1.0033 

1.040 

1 

0.0081 

q 

0(q) 

0(q) 

Upper  b 

Lower  b 

lie  ||% 

0 

1.0742 

1.025 

1.626 

1 

64.51 

l 

1.0482 

1.017 

1.306 

1 

6.21 

2 

1.0203 

1.013 

1.  153 

1 

0. 178 

Realize  that  9(q.  , )/0(q. )  »  -  .  This  is  the  effect  of  the  spectral 

i+1  i  mi+1 

accuracy  of  the  error  estimator.  The  same  effect  is  seen  also  in  the 

4 

improvement  of  the  upper  bound.  This  also  shows  that  the  term  0(d  ) 
neglected  in  (6.2)  decreases  as  q  Increases. 

Let  us  underline  that  in  the  described  example  the  thickness  is  no 
small  (in  fact  when  d  =  2  then  fl  is  a  square)  and  there  is  now  boundary 
layer  in  the  solution. 

7.  The  problem  in  the  presence  of  a  boundary  layer 

Let  us  consider  the  case  f  =  1  where  the  solution  has  a  boundary  layer 

Hence  f  belongs  to  the  class  T„  with  a  large  6.  Table  7. 1  shows  the 

P 

effectivity  index  and  the  relative  error  for  different  q. 
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Table  7.1  The  effectivity  index  and  the  error  for  the 
solution  with  the  boundary  layer 


d  d  =  0.2  .  d  =  1  d  =  2 


a 

a  = 

U.  6 

a  = 

=  1 

a  - 

-  <1 

q 

0 

lie  ||% 

e 

He  ||% 

0 

He  ||% 

0 

1.011 

9.88 

1.063 

42.53 

1.  140 

87.65 

1 

1.46 

0.505 

1.463 

5.  14 

1.454 

12.291 

2 

1.778 

0.  197 

1.778 

2.01 

1.778 

4.736 

3 

1.933 

0.  112 

1.933 

1.  14 

1.933 

2.681 

4 

1.911 

0.079 

1.908 

0.803 

1.912 

1.883 

5 

1.760 

0.064 

1.743 

0.655 

1.760 

1.529 

6 

1.554 

0.056 

1.511 

0.587 

1.554 

1.358 

7 

1.348 

0. 053 

1.292 

0.557 

1.348 

1.269 

We  see  that  the  effectivity  index  is  growing  with  the  value  q,  achieves 
maximum  at  q  »  3  and  then  decreases.  It  is  also  interesting  that  the 
effectivity  Index  is  very  Insensitive  to  the  thickness  d.  The  behavior  of 
0  can  be  explained  as  a  combination  of  the  influence  of  large  /3  and  the 
spectral  accuracy. 

Let  us  now  divide  w  into  subdomains.  In  our  case  u  =  (-1,1)  and  we 
will  divide  u>  into  subdomains  1^  =  (Xj^.x^),  i  =  1....6  and  will  assume 
<p  -  <P^  =  1  on  1^,  =  0  on  Ij,  i  *  J.  This  of  course  is  not  exactly 

in  our  framework  of  assumptions.  Nevertheless  we  can  understand  this 
selection  of  <p  as  an  approximation.  Hence  we  will  speak  about  indicators 
instead  of  estimators.  Let  us  consider  the  6  domains  defined  by  the 
sequence  {x^>  =  {-1,  -0.75,  -0.50,  0.0,  0.50,  0.75,  1.00}.  The  table  7.2 
shows  the  indicators  in  these  domains  with  uniform  q.  Because  of  symmetry  we 
report  the  indicators  on  first  three  subdomains  only. 
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Table  7.2.  The  error  indicators  for  6  domains 


d  = 

=  0.2 

wi 

U3 

II 

q  =  2 

1 

0.576  -1 

0.274  -1 

2 

0.359  -5 

0. 126  -5 

3 

0.224  -9 

0.693  -9 

We  clearly  see  that  the  boundary  layer  is  well  indicated. 

In  Table  7.3  we  show  the  value  of  the  estimator  in  the  first 
subdomain  1^  an  function  of  q 


Table  7.3.  The  error  indicators  In  the  subdomain  1^ 


q 

1 

2 

3 

4 

5 

6 

7 

§(q) 

0.5761-1 

0.2741-1 

0. 1687-1 

0. 1173-1 

0.8768-2 

0.6875-2 

0.5578-2 

8 

9 

0.4644-2 

0.  3945-1 

Figure  5.1  shows  the  graph  of  the  function  £(q).  We  see  g(q)  «  q  p 

1  <  (3  <  2.  Let  us  mention  that  this  rate  is  due  to  the  combined  effect  of 

2 

singularity  in  the  corner  (of  order  r  tg  r)  and  the  boundary  layer. 

Table  7.2  indicates  that  it  is  advantageous  to  use  nonuniform  q.  In 
Tables  7.4a  and  b  we  show  the  error  for  different  distributions  of  q  for 
sets  of  18  and  20  subdomains.  Because  of  symmetry  we  only  list  q  in  the 
first  9  (10)  domains. 

We  used  the  following  two  sets  of  subdomains  defined  by  the  meshes 
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a)  18  subdomains. 


-1.000,  -0.990,  -0.975  -0.950  -0.925  -0.900  -0.850  -0.750 

-0.500  0.00, 

b)  20  subdomains 

-1.000,  -0.998  -0.990  -0.975  -0.950  -0.925  -0.900  -0.850 

-0.750  -0.500  0.00 

The  thickness  d  =  0.2  and  f  =  1.  We  report  in  the  Table  7.4a  the  error, 
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the  effectivity  index  and  W  =  ^  (qg  +  1)  which  (crudely)  estimates  the 


e=l 


work. 


Table  7.4a.  The  distribution  of  q  for  18  subdomains 
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The  Table  7.4b  shows  analogous  results  for  the  case  of  20  subdomains. 
The  first  two  subdomains  are  creating  together  the  first  subdomain  in  the 
case  of  18  elements.  This  is  indicated  by  a  dashed  line  in  the  Table  7.4b. 

Table  7.4b.  The  distributions  of  q  for  20  domains 


case 

Wi 

1 

2 

3 

4 

5 

6 

7 

8 

1 

2 

3 

4 

5 

6 

6 

7 

7 

2 

2 

3 

4 

5 

5 

6 

5 

6 

3 

m 

3 

3 

3 

3 

3 

3 

3 

4 

B 

2 

2 

2 

2 

2 

2 

2 

5 

B 

’  2 

2 

2 

2 

2 

2 

2 

6 

B 

2 

2 

2 

2 

2 

2 

2 

7 

i 

1 

1 

1 

1 

1 

1 

1 

8 

i 

1 

1 

1 

1 

1 

1 

1 

9 

i 

1 

1 

1 

^B 

1 

1 

1 

10 

i 

1 

1 

1 

1 

1 

1 

Error  % 

0.201 

0.  126 

0.  106 

0.094 

0.093 

0.089 

0.92 

0.088 

0 

1.75 

1.74 

1.52 

1.31 

1.21 

1.  13 

1.  11 

1.04 

U 

23 

25 

26 

27 

27.5 

28 

28 

28.5 

<  (q. +  1  )  +  (q, +1 )  y'' 

In  the  Table  7.4b  we  define  ^  W  =  -  +  }  (q^  +  1)  to  make  a 

2  1=1 

comparison  possible.  If  the  degrees  in  and  are  the  same  then  we 

obtain  identical  results  as  in  Table  7.4a. 

We“  see  from  the  Tables  7.4a,b  that  the  optimal  q  distribution  depends 
on  the  error  we  wish  to  achieve  and  on  the  subdomains.  The  main  error  Is  in 
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the  area  of  the  boundary  layer.  Comparing  the  cases  6  and  7  in  Table  7.4b 
we  see  that  the  case  6  leads  to  a  smaller  error  than  case  7.  This  underlines 


the  importance  of  the  right  q  distribution. 

We  have  shown  in  previous  sections  that  the  error  estimator  is  not 
influenced  too  much  by  errors  on  the  subdomains  which  are  not  too  close  to 
the  given  one  (in  the  scale  of  d). 

To  show  the  sensitivity  of  the  indicators  with  respect  to  the 
q-distribution,  let  us  show  in  Table  7.5  the  indicators  in  the  subdomains 

Wg  (in  the  case  of  18  subdomains)  for  some  cases  listed  in  Table  7.4a. 

Table  7.5.  Error  estimators  for  the  q  distribution 

given  in  Table  7.4 


Case 

*3 

*4 

g5 

1 

0. 104-2 

0. 782-3 

0.867-3 

2 

0. 162-2 

0.719-3 

0.  866-3 

4 

0. 163-2 

0.719-3 

0. 866-3 

6 

0. 163-2 

0.719-3 

0. 866-3 

9 

0. 105-2 

0.871-3 

0. 487-3 

Estimators  -  &n  were  completely  Independent  of  the  q-distribution. 
b  y 


8.  Problem  of  adaptive  modelling 

The  goal  of  adaptive  modelling  is  to  design  the  distribution  of  qj  so 
that  for  given  error  tolerance  the  work  is  minimal.  The  quantitative 
characterization  of  the  work  could  be  different.  Here  we  will  use  the  crude 
measure  we  introduced  in  the  Section  7,  namely  W  =  ^ (q^  ♦  1). 

We  will  show  one  of  the  many  possible  adaptive  strategies  in  two  examples. 
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Example  1.  Let  us  consider  the  problem  we  discussed  in  the  previous  section, 

namely  f  =  1  and  d  =  0.2  We  will  use  the  mesh  which  divides  the  domain 

into  18  subdomains.  The  target  relative  accuracy  is  x  =  0.  1%. 

First  we  solve  the  problem  with  q  =  1  (uniform)  and  compute  the 

error  indicators  on  the  particular  domains.  We  single  out  the  subdomains 

where  the  error  is  not  necessarily  made  smaller.  The  principle  is  to  add  Ihe 

O  T2 

smallest  indicators  £7  so  that  their  error  is  =—  B(u,u).  We  conclude  that 

i  b 

q  has  to  be  increased  on  8  subdomains,  namely  and  -  w^g.  The 

goal  is  to  increase  q^  on  these  domains  so  that  we  obtain  the  error  in  the 
prescribed  tolerance.  Using  the  fact  that  the  indicator  on  a  domain  is 

not  essentially  influenced  by  the  indicators  on  the  other  subdomains,  we  will 
"freeze"  the  solution  on  all  subdomains  except  the  one  under  consideration. 

By  this  procedure  we  obtain  the  following  distribution  of  q: 

9, 6, 4, 2, 1, 1, 1, 1, 1  (Because  of  the  symmetry  we  report  here  and  in  what  follows 
only  the  subdomains  for  x  <  0. )  By  this  distribution  we  achieved  the 
relative  error  0.077%.  Comparing  this  distribution  with  Table  7.4a  we  see 
that  the  used  strategy  yielded  the  right  distribution. 

Example  2.  Here  we  assume  that 


and  that  d  =  0.2,  a  *  0. 1. 

We  will  use  10  subdomains  characterized  by  the  mesh  x^  =  -1,  -0.8, 

-0.6,  -0.4,  -0.2,  -0.0,  0.2,  0.4,  0.6,  0.8,  1.0,  with  the  target  accuracy 
0.1%.  In  this  case  the  main  error  is  in  the  place  where  the  solution  is 
unsmooth,  i.e.  in  the  middle  (x  ~  0)  of  the  domain.  The  adaptive  procedure 
gave  the  following  distribution  of  q^  :  1,  1,  1,  2,  3.  The  relative  energy 
error  achieved  was  0.0743%. 
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